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Abstract. We report on an analysis of the polarization of single pulses of PSR B0329+54 at 328 MHz. We find that the 
distribution of polarization orientations in the central component diverges strongly from the standard picture of orthogonal 
polarization modes (OPMs), making a remarkable partial annulus on the Poincare sphere. A second, tightly clustered region 
of density appears in the opposite hemisphere, at a point antipodal to the centre of the annulus. We argue that this can be 
understood in terms of birefringent alterations in the relative phase of two elliptically polarized propagation modes in the 
pulsar magnetosphere (i.e. generalised Faraday rotation). The ellipticity of the modes implies a significant charge density in the 
plasma, while the presence of both senses of circular polarization, and the fact that only one mode shows the effect, supports the 
view that refracted ordinary-mode rays are involved in the production of the annulus. At other pulse longitudes the polarization 
(including the circular component) is broadly consistent with an origin in elliptical OPMs, shown here quantitatively for the first 
time, however considerable non-orthogonal contributions serve to broaden the orientation distribution in an isotropic manner. 
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1. Introduction 

Whenever sufficient sensitivity is available, radio pulsar emis- 
sions are seen to be rich in phenomenology. Their polar- 
ization is no exception to this rule. The dependence of lin- 
ear polarization position angle on pulse longitude (i.e. ro- 
tational phase) can, for some pulsars, be explained as aris- 
ing in the vicinity of the magnetic pole, polarized linearly 
at the position angle of the sky proj ection of magnetic field 
lines (Radhakrish nan & Cookelll969l) . For other pulsars this 
is not the case, and for some of these the distribution of po- 
sition angles (PAs) in individual pulses has been shown to 
be bimodal about two values separated by 90° — so-ca lled 
orthogonal polarization modes (OPMs; e.g . Manchester e^jlJ 
1975; Backer et al. 1976; Cordes et al. 1978; Backer & Rankit 



1980; Stinebri ng et al.lll984aibl) . Many pulsars also show 90° 
"jumps" in their position angle profiles, a fact that received ex- 
planation in the discovery of OPMs through a longitude depen- 
dence of the relative intensities of the modes, which themselves 
tend to have positio n angle swings consistent wi th the magnetic 
pole model of Radhakrishn an & Cookd ([1969). These studies 
of OPMs also found evidence for deviations from orthogonal- 
ity in the fact that PA distributions were broader than expected 
and/or their peaks were not separated by 90°. 

Attempts have been made to explain these deviations from 
orthogonality by means of superposition of two modes that 
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are not orthogonal due to origins on different field lines and 
subsequent birefrin gent refraction (e.g. lStinebring et al.ll9 84a: 
lMcKinnonll200 3a). by means of the superposition of a range 
of modal orientations arising from a distribution of field lines 
that are visible due to their finite beam width ( iGil & Lvna 
1995), or due to the presence of two instantaneously orthog- 
onal modes, the orientation of which v aries with time due 
to coherent wave coupling effects (e.g. ICheng & Rude rman 
Il979t iLvubarskii & Petroval Il999h . Of these, on ly the first 
and thir d scenarios have been quantitatively tested. McKinnon 
(2003a) showed that the distribution of PAs, and the shape 
of the average PA curve of PSR B20 16+28 are consistent 
with the superpositio n of t wo non-orthogonal modes. On the 
other hand, IPetroval J2003I) showed that the PA and circular 
polarization of the average pulse profiles of PSR B0355+54 
and PSR B0 628-28 are co nsistent with the predictions of 
ILvubarskii & Petroval dl999l) for alterations in the polariza- 
tion due to coherent effects. However, neither work examines 
whether the data are also consistent with other models, and in- 
deed the statistics considered in each case are not well suited 
to distinguishing between models. 

The question of the origin of circular polarization in pul- 
sars is also not addressed by the magnetic pole model, in which 
the polarization is expected to be linear due to the high mag- 
netic field strength. The fact that orthogonal linear polariza- 
tion states have been found to be associated with opposite 
signs of circular polarization leads naturally to the sugges- 
tion that the modes are in fact elliptically polarized orthog- 
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onal modes llCordesetalJll978l) . although to our knowledge, 
to date no observational tests of the expected proportionality 
between linear and circular polarization under this hypothesis 
have been made. In explaining the origin of elliptical OPMs, 
most authors point to the so-called polarization-limiting region 
(PLR), where birefringent propagation effects no longer signif- 
icantly alter the relative phase of the modes, as the determi - 
nant of the observed polarization (Cheng & Ruderman 1979). 
It has been proposed that either the propagation modes them- 
selves are elliptica l at this point, requiring a net charge den- 
sity in the plasma i|GiejigjfeRud^ marll979l:lAllen & Melrosd 
1982; von Hoensbroech et al. 1998), or that weak birefringence 
in the vicinity of the PLR causes initially linearly polarized 
rays to suffer changes in their polarization if their position an- 
gle deviates from that of the local linear modes, for exam ple 
due to rotational aberration or refraction (|Chen ^ & Rudermanl 
19791 Lvubarskii & Petrovall999tlPetrova & Lvubarskrlz OOO: 
Petroval200lll2003l) . 

In this work we make a detailed study of the polarization of 
single pulses from PSR B0329+54, employing new techniques 
in the hope of placing greater constraints on the nature and ori- 
gin of ellipticity and non-orthogonality of pulsar polarization. 

2. Methods of Analysis 

2.1. Stokes Parameters of OPM 

Since the focus of this study is the inconsistency of observa- 
tions with a simple model of superposed OPMs, it is pertinent 
to begin with a clear picture of the features expected under 
such a model, before discussing ways in which the observa- 
tions may deviate from it. The distribution of Stokes parameters 
expected under the incoherent superposition of two elliptical 
OPMs, and techniques for reconstructi on of the modal intensi- 
ties h ave been considered in detail by McKinnon & Stinebring 
(2000) .We briefly re-iterate in a more compact vector form be- 
fore considering the fluctuation statistics. 

The modes are taken to be completely polarized, with fixed 
position angles and degrees of circular polarization at a given 
pulse longitude. This means that the vector pi = {QU V) for a 
given mode i always has the same orientation, regardless of its 
length (i.e. the intensity of the mode). The condition of orthog- 
onality of the electric field vectors requires that the polarization 
states be antiparallel in the Poincare sphere (p-space). Since the 
Stokes 4-vector of the incoherent sum of two polarization states 
is simply the sum of their respective Stokes vectors, all resul- 
tant states must lie on the line defined by the modal orientation. 
The observed Stokes parameters at any instant can be written 
as: 

l = h+h (1) 

P = Pi + P2 

= (/i - /a)l», (2) 

where the contributions from the two modes are denoted with 
numerical subscripts, with l m being a unit vector in the direc- 
tion associated with the mode denoted "1". 



The procedure of "mode separation" (determination of 



modal intensities) follows directly as: 

h = d+\p\)/2 (3) 
h = d-\p\)/2, (4) 

while the modal orientation is simply given by 1,„ = p/\p\. Here 



we have explicitly identified mode 1 as the greater in inten- 
sity. Note that \p\ is subject to bias in the presence of measure- 
ment noise. Although this can be corrected, if some alternative 
means is available for determining the modal orientation, since 
components of p orthogonal to 1,„ contribute only noise to \p\ 
one may eliminate the bias and improve sensitivity by substitut- 
ing p \,„ for \p\. Prior information on 1,„ could, for example, be 
obtained from studies of the statistics of p (Sect. 12.3. U . or for 
mode-separation of individual pulses, from the average profile. 

By measuring appropriate statistics of the Stokes vectors, 
information is available on the fluctuation statistics of the 
OPMs. Characterising the latter by the modal variances cr\ and 
cr? and their covariance cr 12 , one observes that 

crj = o\ + cr\ + 2cr 12 (5) 

°"m = °\ + °2 - 2o"i2, (6) 

where cr^ is the variance of p ■ l m . From this the covariance of 
the modes can be calculated as 

o-u = (oj - cr 2 m )/4. (7) 

In practice the measured variances will be biased by the pres- 
ence of measurement noise, which can be estimated using the 
off-pulse covariances (and the pulsar intensity, if it contributes 
appreciably to the system temperature), and subtracted. 

2.2. Characterising the Distribution on the Poincare 
Sphere 

The defining feature of elliptical OPMs is the predicted con- 
stant orientation of p, regardless of pulse-to-pulse fluctuations 
in / and \p\. This fact is put to good use in the display of 
histograms of the position angle of linear polarization, iff = 
i tan -1 U/Q, which at a given longitude can take one of two 
allowed values offset by 90°. In the presence of instrumental 
noise these broaden to give the characteristic bimodal distribu- 
tion of OPM. The position angle distribution gives information 
on the orthogonality of the linear component of the polariza- 
tion, but the question must also be asked whether the circular 
component is as expected under elliptical OPM. Attempts to 
answer this question have been hampered in the p ast due to 
the use of inappropriate statistics. McKinnon (2002) notes that 
the picture of elliptical OPM is consistent with observations 
showing distributi ons of Stokes V/I that are broad and cen- 
tred near zero jBacker & Rankinll980HStinebring et al.ll9 84a) 
(or eq uivalently, we note, that (|V|) tends to greatly exceed 
| (V) |; Karastergiou et al. I l2003l) . However, such a distribution 
could not be considered a particularly distinctive feature of 
OPM and could easily be produced even by a mechanism for 
the production of circular polarization that is independent of 
the linear OPM phenomenon. A much more stringent test is 
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to check that the observed V is, along with Q and U, con- 
sistent with a constant orientation of p. The natural comple- 
ment to the position angle in this regard is the ellipticity angle, 
X — \ tan -1 V/ ~<JQ 2 + U 2 : together they completely specify the 
orientation of p via the spherical coordinate angles 2tfr, 2x- The 
distributions in both parameters should be bimodal under OPM 
(unless one mode always dominates), unlike the distribution of 
V/I, which may be unimodal. 

While a measurement of the joint probability density func- 
tion, /(iA,^), contains sufficient information to detect the pres- 
ence of elliptical OPM, this choice of parameterisation is not 
ideal. This is because a given solid angle element on the 
Poincare sphere subtends an area in (^,^-space that itself de- 
pends ovix- Specifically, 



da = 4di^d;f cos 2^. 



(8) 



This means that points closer to the equator (x — 0) are given 
greater weight in the distribu tion than j ustified by their density 
on the sphere (see also McKinnon 2003b). This problem could 
be circumvented by using f(i[f,x)/ cos 2y instead, but the prob- 
lem remains that the spatial scale with which features in the 
distribution are represented in Cartesian ijj,x coordinates (for 
example in a plot) varies with x- 

A solution to this problem is to measure f{2\]/, sin 2x), 
which satisfies the equal-area condition 



da = d(2i/r)d(sin2^). 



(9) 



This choice of parameterisation is known as Lambert's cylin- 
drical equal-area projection. In fact, any equal-area projection 
has the desired properties of providing true solid-angle densi- 
ties and representative sizes for features. Lambert's cylindrical 
equal-area projection suffers from severe shape distortions at 
both poles (Q = U = 0, i.e. x = ±45°) and we shall only make 
use of it in the case where one of the two dimensions is to be 
averaged in order to use the other along with pulse longitude 
in a two-dimensional display. That is, we supplement the tra- 
ditional iff versus longitude density display with a sin 2x versus 
longitude display that gives information on the distribution of 
ellipticity as a function of longitude. 

To examine the distribution of polarization orientations at a 
single pulse longitude, it is desirable to have minimum (zero) 
distortion at both of the orientations associated with the OPMs. 
A cylindrical equal-area projection with its equator containing 
the OPM orientations would suffice, however the nature of de- 
viations from OPM in PSR B0329+54 we describe below mo- 
tivates the use of a projection where the distortion is axisym- 
metric about the modal orientation. The only such projection to 
also have the equal area property is Lambert's azimuthal equal 
area projection, f(x, y) where 



x = pcosA 
y = psinA 



(10) 
(ID 
(12) 



and A, 8 are the azimuth and latitude of p in a rotated frame with 
l m as its "north" pole. This projection maps the entire sphere to 
a disc, with the "north" pole in the centre and the "south" pole 



distorted to a ring about the circumference. In order to avoid 
this severe distortion, we truncate the projection at 8 < and 
display the southern hemisphere in an adjacent projection with 
-l m as its central pole. 

2.3. Testing for Orthogonality 

2.3.1 . Eigenanalysis of the Covariance Matrix 

Previous studies of pulsar polarization have frequently revealed 
position angle distributions with components that are signif- 
icantly br oader than that e xpected due to instrumental noise 
alone (e.g. Stinebrin g et alJll984ah . indicating the presence of 
some non-orthogonal radiation. Deriving as it does from the 
ratio of two noisy quantities (Q and U), the linear polariza- 
tion position angle possesses complicated statistics, depending 
not only on the statistics of measurement noise but also on the 
statistics of variations in the linear polarization of the source. 
Although much ef fort has been expended modelling them (e.g. 
lMcKinnonl200 3b). it is clear that as a means of detecting non- 
orthogonality, the PA distribution is less than ideal. The distri- 
bution of ellipticity angles suffers from the same problems, in 
addition to those discussed in Sect. I2.2I We suggest instead a 
more robust test of the consistency of polarization observations 
with an origin in OPM. 

From the fact that OPMs contribute to the observed values 
of p only along a single vector, it follows that a simple and di- 
rect means of detecting non-orthogonal emission is to check the 
statistics of the components of p perpendicular to l m . In order 
to do this, it is necessary to define a new orthogonal basis for 
p, which has 1,„ as one of the basis vectors. The question then 
arises of what to use for a working value of l m , since when the 
orthogonality of all states is in question, the mean polarization 
vector may not be parallel to l m . Under the assumption that the 
majority of the fluctuations in p are directed along l m , a sen- 
sible choice for this vector is the direction of greatest variance 
in p, which is also the least-squares estimate of the direction 
of fluctuat ions. The met hod of Principal Components Analysis 
(PCA; e.g. lJollifell 986) is suggested as a means of finding this 
vector. PCA is based on the fact that the set of eigenvectors of 
the covariance matrix of a multivariate statistic represent an or- 
thogonal basis in which the variations in the different axes have 
zero covariance 1 . That is, given the covariance matrix 



C = ((p-{p))(p-(p)) T ) 



(13) 



where p is a column vector, there exists a set of 3 orthonormal 
vectors e,- (eigenvectors) and scalars A/ (eigenvalues) such that 



Cei = A-fii. 



(14) 



It is easily shown that the covariance matrix is diagonal if the 
vectors are expressed using the eigenvectors as a basis. That is: 



C = {[M(p-(p))][M(p-(p))] T ) 



(15) 



1 We note that a similar tech nique was independently described for 
OPM use by McKinnon (2004) while this paper was in the final stages 
of preparation. 
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A 2 

A 3 



(16) 



where M is a matrix with the eigenvectors as its rows. 

Clearly, the eigenvalues are equal to the variances in the 
components aligned with the corresponding vectors. The eigen- 
vector with the greatest eigenvalue corresponds to our choice 
of 1,„, while the other vectors and their associated eigenvalues 
present a convenient basis for detection and characterisation of 
non-orthogonal radiation. In theory there may be fewer than 
three distinct eigenvalues indicating axisymmetry (two eigen- 
values) or isotropy (one eigenvalue) in the directional variance, 
although the presence of measurement noise makes this have 
zero probability in practice. However, measurement noise also 
introduces bias to the covariance matrix that should be cor- 
rected by subtraction of the covariance of the noise, estimated 
from off-pulse longitudes and including also the contribution 
of the pulsar to the system temperature if significant. Where 
the intrinsic covariance is of a similar magnitude to the un- 
certainty in the off-pulse covariance, this may result in neg- 

1 II 

ative eigenvalues, which presents a problem if A. is to be 
used as a measure of the scale of intensity fluctuations, and 
it must be accepted that estimates will be unavailable in some 
longitude bins. This problem is also familiar from, for exam- 
ple, bias-corrected estimates of the linearly polarized intensity, 
L = (Q 2 + U 2 - a* - a%)^ 2 . 

Having measured the variance of fluctuations in compo- 
nents parallel to and perpendicular to the modal orientation, 
it is useful to define a measure of the degr ee of deviation from 
purely linear fluctuations. After ICloude & Pottierl i 19951) . we 
define the polarization entropy as: 



H=-Y j P i log 3 P i , Pi = 

1=1 



(17) 



where At is the i-th eigenvalue. This quantity is confined to the 
range [0,1], where H — corresponds to fluctuations com- 
pletely confined to one direction, and H = 1 corresponds to an 
isotropic distribution of p - (p). 

2.3.2. Effects of Scintillation 

The methods of the previous section work from the covari- 
ance of the observed signal, which in practice may reflect not 
only the intrinsic variations of the pulsar, but also the effects 
of the interstellar medium on the propagating signal. Of princi- 
ple importance is the variable effective gain of the interstellar 
medium induced by two types of propagation effect, refractive 
and diffractive scintillation. For most pulsars, the time scale 
for variations due to refractive sci ntillation is long (~ 15 days 
for PSR B0329+54: rStinebring etafl feoOO). so the effects can 
be neglected over short observations. Time scales for diffrac- 
tive scintillation, on the other hand, are much shorter. In some 
cases the time scale may be long enough that one can simply 
limit the analysis to a segme nt of time over which the flux den- 
sity is constant (e.g. McKi nnonl2004 . but in other cases such 
intervals do not last long enough to sensitively obtain repre- 
sentative statistics. With a diffractive scintillation time scale of 



~ 148 s at 328 MHz dCordesl H986). this is the case for the 
observations of PSR B0329+54 reported here. 

In order to separate the (co-)variance induced by scintil- 
lation from that int rinsic to the pulsar, we emp loy fluctuation 
spectral techniques (Edwards & Stappers 2003). We define the 
(longitude-resolved) polarization spectral density tensor as the 
element-wise product of the (discrete, vector) Fourier trans- 
form of pj - (j>j) with its Hermitian transpose, where / indexes 
pulse number and the longitude dependence is implicit as else- 
where in this work. That is, 



S k = PkPk , where 

N 



(18) 
(19) 



N is the number of pulses, and 1 denotes the 
Hermitian transpose. Here £ is a frequency index (correspond- 
ing to a frequency of k/N cycles per pulse period), which is 
most easily interpreted if it is taken in the range (— N/2, N/2], 
For convenience, in practice the observation can be broken into 
segments for which S is computed, and the results averaged. 
The tensor S contains a wealth of information about the polar- 
ization fluctuations, offering among other possibilities, a means 
for distinguishing different sources of variation in p via their 
different fluctuation frequency structure, and a chance to form 
the complex covariance matrix of the analytic signal and check 
for periodic patterns of elliptical (rather than only linear) form 
(Edwards 2004, in preparation). In the present context we sim- 
ply make use of the fact that the sum of St over all k gives, via 
Parseval's theorem, the covariance matrix (Eq. 1131 . while the 
covariance matrix of a bandpass-filtered version of the data can 
be computed as a sum over a suitably limited range of k. 

The effect of scintillation on the observed signal pj is a 
time-varying, multiplicative gain. The observed signal can be 
written 



pj = (l+mpjjp*., 



(20) 



where m is the modulation index due to scintillation, pj is a 
realisation of a random process of zero mean and unit variance, 
and p l is the intrinsic signal. Breaking the intrinsic signal into 
a mean contribution plus a variable contribution, Eq.[20]can be 
expanded to give 



Pj -</>> = (l + m Pj ) (p*. - (p'yj + mpj (p f ) . 



(21) 



Assuming that the instrinsic and scintillative variations are in- 
dependent, the covariance matrix (Ea.ll3> becomes 



C = (l+m 2 )C +m 2 (p i )(p i ) T . 



(22) 



This shows that the uncorrected covariance matrix is biased and 
scaled by the presence of scintillation 2 . Replacing multiplica- 
tions in the time domain of Eq. with convolutions in the 



2 As an aside we note that the presence of scaling implies, surpris- 
ingly, that scintillation can improve sensitivity to intrinsic variance, 
even when the data contain an equal balance of scintillative amplifica- 
tion and deamplification events. 



R. T. Edwards and B. W. Stappers: Ellipticity and Deviations from OPM in PSR B0329+54 



5 



frequency domain, one finds that the vector spectrum of Eq. 
[^consists of the sum of the spectrum of intrinsic variations 
from the mean, its convolution with the spectrum of the scin- 
tillation variations, and the convolution of the spectrum of the 
scintillation variations with the spectrum of {^), i-e. a scaled 
sine function at zero frequency. The condition of independence 
implies that this sum carries over to Ea. ll8li ust as it did to Eq. 
l22l The spectral response of the scintillation variations has a 
characteristic width 1 /t (where r is the scintillation timescale), 
corresponding to Ak = NP/t, where P is the pulse period. In 
the usual case that P/t «: 1, we see that the bias term in the 
observed covariance is confined to a small fraction of the spec- 
tral density tensor around zero frequency. Under the reasonable 
assumption that the fraction of intrinsic fluctuation energy at 
these frequencies is small, this region can simply be omitted in 
the sum over k that is performed to form the covariance matrix, 
to eliminate the bias. The remaining scale factor 1 + m 2 can 
be obtained (given (p)) from the complementary sum of the 
spectral tensor over small \k\, which is over-constrained, ne- 
cessitating a fit. In practice, since m and (to first order) pj are 
independent of pulse longitude, it can be obtained easily and 
accurately by adding all on-pulse longitude bins in the Stokes 
I signal and extracting the analogous quantity from the fluctu- 
ation power spectrum of the result. 

3. Observations and Discussion 

3.1. Observations 

We used the Westerbork Synthe sis Radio T elescope with the 
PuMa pulsar backend ( Vo ute et alJ I2002I) to observe PSR 
B0329+54 in a 10 MHz band centred at 328 MHz. PuMa was 
configured as a digital filterbank, producing samples in all four 
Stokes parameters over 128 frequency channels with a sam- 
ple interval of 409.6 /us. We found that our standard digitisa- 
tion using 2 bits per sample per channel was insufficient to 
avoid clipping the brightest pulses, and caution that this may 
also have been a problem in previously published polarime- 
try of this and other bright pulsars. This problem was avoided 
by re-observation of the pulsar with 8 bits of digitisation. In 
offline analysis the data were corrected for instrumental po- 
larization effects determined using the procedure described in 
the Appendix, followed by removal of the frequency-dependent 
position angle rotation caused by interstellar Faraday rotation, 
and samples were summed across all frequency channels after 
correcting for delays due to interstellar dispersion. The resul- 
tant time series was divided into segments corresponding to the 
apparent pulse period to give an array of 16400 pulses in 1744 
pulse longitude bins, 220 of which were used in further analy- 
sis. 

This pulsar is known to exhibit at least four mai n modes o f 
emission, each with a different pulse profile ( Ba rtel et alJl9 82 ), 
two of which are commonly seen at low frequencies. By form- 
ing pulse profiles in sub-integrations of 100 pulses, we deter- 
mined that the pulsar underwent a mode change at pulse num- 
ber ~ 12800. The profiles formed by adding pulses 0-12500 
and 13000-16400 are consistent with the so-called "normal" 
and "abnormal" modes reported by previous authors. 



3.2. Distribution of Polarization Orientations 

As a first step in the characterisation of the distribution of ori- 
entations of the polarization vector p, we formed frequency of 
occurrence histograms in position angle and sin 2x (Sect. 12.21 
as a function of pulse longitude. To reduce the effect of noise, 
only those samples with \p\ > 10cr p , where cr p is the quadra- 
ture mean of the variance measured in Stokes Q, U and V in 
off-pulse longitude bins, were included in the histograms. The 
result for the normal profile mode is shown in Fig.^ 

The histograms for the abnormal mode were qualitatively 
very similar to those of the normal mode, upon which the 
remainder of the discussion will focus. However, we note 
as an aside that the position angle distributions appeared 
offset between the two profile modes. Although the pres- 
ence of non-OPM means that the magnetic pole model of 
Radhakrishna n & Cookd i 19691) is not stri ctly applicable, for 
the purpose of quantifying the offset, after lOil & Lvnel £l995) 
we attempted to fit the model to position angles determined us- 
ing local maxima in the histograms near to a curve made by eye 
to approximately fit the less distorted of the polarization modes 
(Fig-D- We found that for the abnormal mode, a fit to the fidu- 
cial position angle parameter 1^0 while using values determined 
from the normal mode for other parameters performed much 
better than a fit where if/Q was fixed and all others allowed to 
vary (RMS deviation 3.13° versus 3.91°), and was comparable 
to the best fit with all parameters allowed to vary (RMS 3.06°). 
We therefore conclude that the offset is consistent with a shift 
in position angle of ~ -5.5°, rather than a change in the appar- 
ent viewing geometry or a longitude offset as might be induced 
by differential aberration and retardation. 

While the position angle distribution of Fig.^is consistent 
with previously publishe d results of lower resolution showing 
quasi-orthogonal modes (Gil & Lvne 1995), the ellipticity dis- 
tribution shows features of a kind never seen before in any pul- 
sar, owing most likely to the fact that previous studies have 
used V/I instead of ellipticity, causing OPM-related features to 
be washed out due to fluctuations in \p\/I. Most striking is the 
strong right-circular polarization seen under the main central 
component, which has no corresponding left-circular compo- 
nent of equal ellipticity as would be expected if the circular 
polarization is due to the OPM clearly seen in the position an- 
gle distribution (Sects. |2~T1 & 12. 2> . Also of interest is that the 
trailing component (longitude ~ 64°) has a distribution that 
is roughly bimodal about the zero line, as expected under el- 
liptical OPM, while the leading component (longitude ~ 41°) 
appears to have a unimodal ellipticity distribution. Also appar- 
ent is that polarized emission is occasionally detected in the 
vicinity of pulse longitudes 33° and 70°, corresponding to the 
additiona l emission compon ents detected in total intensity by 
iGangadhara & GuptalfeOOll) . 

Much more intriguing behaviour is made apparent when the 
full two-dimensional orientation distribution is considered for 
particular longitude ranges. In Fig.|5]we display these distribu- 
tions averaged over several longitude intervals, using the pro- 
jection described in Sect. 12.21 In what follows we refer to the 
modes occurring in the left and right halves of each projection 
as modes "1" and "2" respectively. Addressing the distributions 
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Fig. 1. Longitude-dependent polarization behaviour of PSR B0329+54 in its normal profile mode at 328 MHz. Plotted are the 
mean total and polarized intensity (black and red, top panel) and histograms of position angle (middle panel, plotted twice for 
continuity) and sin 2x where x is the ellipticity angle, along with the corresponding parameters of the mean polarization vector 
(black lines). For convenience the ordinate scale of the bottom panel is labelled non-linearly in terms of x- The histograms are 
normalised by the peak density and plotted in the colour scale shown at the bo ttom. The dashed line in the middle panel shows 
the prediction under the magnetic pole model ( Radha krishnan & Cookelfl969i) . see text. The smooth behaviour seen is intrinsic 
to the pulsar, the longitude bin spacing is determined by the instrumental resolution. 
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Pulse longitude (") 

Fig. 2. Distribution of polarization orientations in eight longitude intervals, marked in the pulse profile at the bottom and shown 
in left-to-right, top-to-bottom order. Each panel shows the distribution in Lambert's azimuthal equal-area projection with poles 
set near to the typical orientation of the mode showing the least scatter. The projection is interrupted at the equator and plotted in 
two hemispheres, to give low distortion near both poles, marked by asterisks. Lines of constant iff andx are marked for reference; 
meridians (constant if/) converge atx = 45° (pure right-circular polarization) near the top of the projection. The colour scale is 
shown in Fig.^ eacn density being normalised by the peak in each longitude interval. 
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in longitude order, we see that the leading component is con- 
sistent with purely linear OPM in mode 2, while by pulse lon- 
gitude 50° the modes have switched in dominance and become 
somewhat elliptical. As the pulsar rotates, mode 2 begins to in- 
crease again in strength, and apparently has a greater spread 
in its orientations than mode 1. Over the course of the cen- 
tral component the distribution associated with mode 2 deforms 
into an arc and eventually an almost complete annulus, while 
mode 1 remains tightly distributed around an elliptical orien- 
tation and eventually concedes dominance to mode 2. Finally, 
in the trailing component the modes are of comparable strength 
and distributed tightly around orthogonal elliptical orientations. 
We discuss our interpretation of this remarkable behaviour in 
Sect. l3.4l but first discuss the remaining observational results. 

3.3. Covariance Analysis 

The results of the previous section leave no doubt that the cen- 
tral component shows strong deviation from the behaviour ex- 
pected under OPM. The case of the leading and trailing com- 
ponents is more difficult to assess due to the fact that the ex- 
pected spread in orientations under instrumental noise depends 
in a complicated way on the distribution of \p\. Instead, we 
used the method of eigenanalysis described in Sect. 12.3.11 In 
our case the dispersed pulsar signal contributes at most about 
one fifth of the total system temperature, justifying the use of 
a single, longitude-independent correction to the spectral den- 
sity tensor for the off-pulse noise, under the caveat that a very 
small (<x =s 0.003 normalized flux units) amount of measure- 
ment noise contaminates the variances for the central compo- 
nent. The characterstic frequency corresponding to scintillation 
on the diffractive time scale (148 s: ICordesl [l986) is — 1 /200 
cycles per period. To ensure the response was eliminated over 
its full frequency extent, we excluded elements of the spectral 
density tensor with \k\/N < 1/50 when computing the covari- 
ance matrix. Using power from \k\/N < 1/50 in the fluctua- 
tion spectrum of the pulse energy, we measure a modulation 
index of 0.16 due to scintillation, in agreement with the mea- 
surements and empirical model of Cordes ( 1986), given our ob- 
serving band. This value was used to correct the overall scale of 
the covariance matrix. The results of the eigenanalysis of this 
matrix are shown in Fig. [3] 

Beginning with the polarization entropy (Eq. I17> . we see 
that the polarization is most disordered in the central compo- 
nent, as one might expect from the distributions seen in the 
previous sections, but still shows detectable entropy in all other 
pulse longitudes accessible to measurement. That the diver- 
gence from pure OPM is significant is confirmed by the fact 
that significant variance is detected in the second and third 
eigenvalues under every component. The fact that the second 
and third eigenvalues are nearly equal in all components ex- 
cept the central peak indicates that, for these pulse longitudes, 
the deviations from OPM show no preferred direction, and 
cannot be caused by position angle distortions or a random 
circular component alone. This is consistent with the analy- 
sis of PSR B1929+10 and PS R B2020+28 at 1404 MHz per- 
formed by McKinnon (2004), who suggests the superposition 




.45 I 1 1 1 1 1 1 i—l 

35 40 45 50 55 60 65 70 
Pulse longitude (*) 

Fig. 3. Results of eigenanalysis. The top panel shows the av- 
erage intensity profile (thick solid line), the square roots of 
the eigenvalues (solid, dashed, dotted thin lines, in descend- 
ing order of value), and the polarization entropy (Eg. 1 171 thick 
dotted line). The middle panel shows the position angle of the 
mean polarization vector (thick line, plotted repeatedly at off- 
sets of 180° for continuity) and of the first eigenvector (thin 
line, plotted repeatedly at offsets of 90° for continuity under the 
equivalence of antiparallel vectors as eigenvectors). The bot- 
tom panel shows the ellipticity angle of the mean polarization 
vector (thick line) along with the ellipticity angle of the first 
eigenvector, plotted twice with opposite signs. 



of (isotropic) randomly polarized radiation (RPR) as the cause. 
On the other hand, in the central component of PSR B0329+54 
all three eigenvalues are significantly different, and indeed the 
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analysis of the directional distribution in the preceding section 
shows that the distribution of p cannot be ellipsoidal. This im- 
plies that the deviations are themselves associated with the pro- 
duction of OPM, as discussed further below. 

It is also interesting to note that there is a suggestion of cor- 
respondence between transitions in the mean modal dominance 
(e.g. Fig.^ an d local minima in the polarization entropy (Fig. 
|3jl. This appears to be the case around pulse longitudes ~ 42°, 
~ 52°, and ~ 63°, however the correspondence is not exact, 
particularly in the leading component. If the trend is in fact 
real and confirmed in other pulsars, it would require that any 
mechanism for the production of two OPMs predicts that under 
conditions leading to OPMs of similar intensity, OPM-related 
fluctuations dominate more strongly over the randomly polar- 
ized fluctuations than elsewhere in the pulse profile. This could 
be the case if the modal intensities tend to be more variable or 
more negatively correlated (Eq.|6j, and/or the randomly polar- 
ized component is weaker or less isotropic. A detailed study 
of a larger sample of pulsars would be necessary to distinguish 
between these possibilities. 

We also note that a smooth position angle curve can be 
constructed from the eigenvector corresponding to the largest 
eigenvalue, in contrast to the position angle of the average po- 
larization vector, which shows gradual 90° transitions rather 
than sharp flips as would be needed for reconstruction of a con- 
tinuous smooth curve. Also the ellipticity angle curve of the 
first eigenvector avoids the problem seen at longitude ~ 39°, 
where near complete cancelling of the linear contributions 
of the OPMs is not accompanied by cancelling of the circu- 
lar component, giving rise to a "spike" in x where (p) / \{p)\ 
sweeps over the left-circular pole (indicating, incidentally, non- 
orthogonal modes, or a consistent, superposed left-circularly 
polarized component). These properties will likely make eige- 
nanalysis a useful technique for detecting polarization fluctu- 
ations driven by OPMs and determining the longitude depen- 
dence of their polarization orientations, even when the signal- 
to-noise ratio is insufficient to detect individual pulses. 

3.4. Interpretation of Strong Deviations from OPM 

The results of the preceding sections show that, while at some 
pulse longitudes the circular polarization is typically in propor- 
tion with the linear polarization such that two clusters are pro- 
duced in the orientation distribution at antipodal points on the 
Poincare sphere, at no pulse longitude is the scatter in orienta- 
tions consistent with instrumental noise alone. Moreover, near 
the peak of the average profile, one mode shows extreme diver- 
gence from the expected orientation, with a correlation between 
position angle and ellipticity that is complex in form. The shape 
is consistent with a broad, incomplete annulus centred upon the 
point opposite to the orientation of the other mode. 

A possible origin for this behaviour lies in birefringent ef- 
fects in the magnetosphere. Specifically, the annular form is 
suggestive of a propagation effect whereby the polarization 
state of incoming rays as represented on the Poincare sphere 
are rotated by a time-varying angle about the axis defined by 
the central point of the annulus. Such an effect is expected if 



the observed radiation passes through a region of plasma where 
the natural propagation modes are different to the ray polariza- 
tion, and a net phase delay occurs between the components of 
the electric field in each of the two modes due to their differ- 
ent group velocities. Thi s effect has been termed G eneralised 
Faraday Rotation (GFR; iKennett & Melrosdfl993) and is fa- 
miliar from, but not theoretically limited to, ordinary Faraday 
rotation in the interstellar medium (about Stokes V, due to the 
circular modes of cold, non-relativistic magnetised plasma), 
and from the effect of retardation plates (rotation of p about 
a linear orientation defined by the optical axis of the material). 
In the case of the pulsar, the polarization of the plasma modes 
can be identified with the slightly elliptical polarization states 
that appear antipodal on the Poincare sphere at the centre of 
the annulus and at the typical orientation of states dominated 
by the other, well-behaved mode. The incoherent superposition 
of radiation in the other polarization mode, which apparently 
does not suffer this effect, would cause the annulus to broaden 
outwards, helping to explain the spread of the observed distri- 
bution. 

This kind of effect was pred icted for pulsar magnetospheres 
bv lChen g & Rudermanl \ 19791) and g iven a quantitative treat- 
ment bv lLvubarskii & Petroval 41999). In their formulation the 
change in polarization is effected in the vicinity of the polariza- 
tion limiting region (PLR; Sect. [3, where the plasma density is 
insufficient to cause total decoherence of the modes, yet densi- 
ties are still high enough to cause significant phase delays be- 
tween the modes. Radiation enters this region as an incoherent 
mixture of the local plasma modes, but due to changes in the 
modal orientation along the ray path caused by the rotation of 
the magnetosphere, each of the incoming rays acquires compo- 
nents in both of the propagation modes, which propagate at dif- 
ferent speeds and alter the polarization of the ray accordingly. 
This p icture deviates from our obs ervations in several ways. 
Firstly, Lvubarskii & Petroval i 19991) assume linearly polarized 
propagation modes, whereas the observations indicate ellipti- 
cal modes, implying plasma with a net charge density rather 
than a pure pair plasma. Secondly, the predicted effect is not as 
simple as rotation about a given axis, for the modal polariza- 
tion orientation, and thus the axis of rotation, varies along the 
ray path. The presence of an annular shape, as expected from 
a near-constant polarization of the propagation modes, may 
therefore place some constraints on the size of the region of 
the magnetosphere contributing significant, variable amounts 
of GFR. Finally, the effect should o nly be capable of indu cing 
one sense of circular p olarization ( Radhakrishnan & Rankin 
1990; Lvubars kii & Petroval 19991) . and should affect both rays 
equa lly apart from a reversal in the sense of circular polariza- 
tion JPetroval200ll) . 

An alternative cause of the misalignment of the polariza- 
tion of the propaga tion modes and the incoming rays, is refrac- 
tion. |P£jr£>xiij^L^uEaS3 i|2QQ3) show that, while the extrao- 
dinary mode propagates under a vacuum dispersion law, the 
ordinary mode can suffer from significant refraction, which, 
under a non-axisymmetric plasma distribution, can cause it to 
move out of the plane of the magnetic field line from which 
it originat ed (and obtained its in itial polarization). The calcu- 
lations of Petrova & Lvubarskii (2000) show that the subse- 
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quent alteration of the polarization state at the PLR can pro- 
duce either sense of circular polarization, as seen in our ob- 
servations. Moreover, since the extraordinary mode is immune 
to refraction, it should not suffer the same PLR effects, con- 
sistent with the tight, centrally peaked distributions of orien- 
tations observed here in mode 1 . Should one o f the modes b e 
produced by conversion from the other, as in IPetroval J200ll) . 
this would imply that refraction occurs above the conversion re- 
gion. An alternative means of producing mode-dependent PLR 
effects is the invocation of [anti-]correlation between the effi- 
ciency of conv ersion and the physical conditions in the PLR 
(Petroval Eool . however for this to be the case the correla- 
tion must be very strong, given the complete absence of an 
annulus in the distribution of states apparently dominated by 
mode 2. That refractio n-driven PLR effects only occur close 
to the magnetic axis dPetrova & Lvubarskiil 120001) is another 
prediction borne out by the observation that only the central 
profile component shows the annular distribution. Many pul- 
sars show strong me an circular pol arization in central, so-called 
"core" components (Rankin 1983 ), which tends to sh ow a cen- 
tral sens e reversal ( Radhakrishnan & Rankin 1990, alt hough 
see also Han et al. 1998) taken bv lPetrova & Lvuba rskii (2000) 
as support for their model of the refraction-driven PLR effect. 
The probable direct detection of this effect in PSR B0329+54 
opens the possibility of good tests of the model through ap- 
plications of the techniques used here on a larger sample of 
pulsars with and without "core" components, and examination 
of the frequency dependence. 

Detailed modelling of this effect is beyond the scope of this 
work, however to prove the basic assertion that GFR can pro- 
duce the spread of orientations observed, we have performed 
some basic numerical simulations. We simulated the observed 
polarization vector as the sum of three components, the ordi- 
nary ray, the extraordinary ray and an RPR component: 



P = Po+Pe+ Pr- 



(23) 



We assume that certain parameters are independently dis- 
tributed according to chosen "reasonable" distributions. These 
are merely intended to provide variation in order to produce 
the qualitative features of the model, and are not expected to 
necessarily closely resemble the distributions under a complete 
physical model. The ordinary mode was presumed to have an 
initially linear polarization orientation with a position angle 
drawn from a von Mises distribution: 



/(2«A) = 



e 



(24) 



(where Iq is a modified Bessel function of the first kind) with 
fj. v = 50°, k = 20, while its (polarized) intensity was drawn 
from an exponential distribution: 



/(/) 



/<0 



(25) 



with ix e = 1 . To model the effect of GFR, the initial state for the 
ordinary mode was rotated by an angle about an axis l m ori- 
ented with ip = 90°, x - 4°, where 9 was also drawn from a von 
Mises distribution (ju v = 90°, k = 1), to give the computed p . 




Fig. 4. Distribution of polarization orientations deriving from a 
simulation involving GFR (see text). The projection used is as 
in Fig. |3 plotted in a linear grey density scale. 

The extraordinary mode (p e ) was presumed to be always ori- 
ented with ift = 0°, x — - 4°, with a (polarized) intensity drawn 
from an exponential distribution of mean fi e = 0.2. The contri- 
bution of RPR (p T ) was calculated by drawing its components 
Q, U, V from independent Gaussian distributions: 



1 



-V2^r 



(26) 



with yu g = and cr = 0.05. This model was realised 10 6 times 
to produce the distribution of orientations depicted in Fig. 0] 
Although clearly the chosen distributions are not perfect, the 
resemblance to the observed distributions is striking, confirm- 
ing the possibility that GFR is in effect. 

4. Conclusions 

We have performed a detailed analysis of the polarization of 
individual pulses from PSR B0329+54 at 328 MHz. The use of 
statistics involving the ellipticity angle instead of Stokes V or 
V/I enabled us to reveal for the first time the rich structure of 
the circular polarization distribution, which is correlated with 
the position angle distribution. 

Displaying the distribution of polarization orientations in 
an equal-area projection, a remarkable structure was uncov- 
ered for pulse longitudes near the centre of the average pulse. 
The annulus-like form of the distribution in one mode is, in 
our view, indicative of Generalised Faraday Rotation (GFR) in 
the pulsar magnetosphere, while the fact that it apparently only 
affects one polarization mode, and is capable of producing ei- 
ther sense of circular polarization, is taken as an indication that 
refraction of the ordinary mode is taking place. The apparent 
ellipticity of the polarization state upon which the annulus cen- 
tres and orthogonal to the typical polarization of states domi- 
nated by the other mode are taken to indicate that the GFR oc- 
curs in a region with elliptically polarized plasma propagation 
modes, indicating a net charge density in the plasma. 

Through the analysis of the covariance of the Stokes param- 
eters, by means of eigenvector decomposition, we have shown 
quantitatively for the first time that the circular polarization 
is consistent with an origin in elliptically polarized orthogo- 
nal polarization modes. Moreover, we find that a significant 
apparently randomly polarized component dilutes the purity 
of OPM-driven variations in polarization state, as also found 
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recently bv iMcKinnonl ( 12004 for PSR B1929+10 and PSR 
B2020+28 at 1404 MHz. 

In this work we have shown how to detect and charac- 
terise deviations from OPM in more powerful ways than previ- 
ously available. Future application of these techniques to other 
pulsars over a broad frequency range should allow renewed 
progress in the resolution of decades-long debate about the ori- 
gin of pulsar polarization. 
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Appendix A: Polarization Calibration 

The Westerbork Synthesis Radio Telescope (WSRT) is an array of 14 
dishes, each of two nominally linear dipoles, which is configured as 
a tied array during pulsar or VLBI observations. Assuming that the 
response of each telescope to the incident electric field is linear, the 
observed voltages can be described using a Jones matrix iJoneJ 19411) 
for each dish, i: 

v, = he. (A.l) 

The procedure for forming a tied array is to add the voltages of all 
telescopes after correcting them for relative phase (due to the source 
direction): 

N 

V,„, = J] e ' WiJ ' e - ( A - 2 > 

i=l 

where A' is the number of telescopes. It is therefore trivially true that 
the tied array signal can be calibrated by determining an overall sys- 
tem Jones matrix, i.e.: 

vtot = /tote, (A.3) 

N 

/.o. = Yj e ~' BiJi - (A - 4) 

The problem of tied array calibration is thus no different to the calibra- 
tion of single-dish observations, for which several recent treatments 
are available that avoid potentially dange r ous sim plifying assumptions 
of ea r lier works illamaker et all|l996l : ISault et alj Il996t lHamakei 
2000; iBrittodExMi Ivan Stratenl 12003 : Ijohnstonl 120021; Ivan Strater 
2004). These approaches are based on the Jones matrix formalism, 
which has seven degrees of freedom in the instrumental model of each 
frequency channel (the eighth, absolute phase, is lost in the compu- 
tation of th e self coherency matrix for single dish or tied array ob- 
servations), [nimaj^ j20(Xj) divides these into a "polrotation", which 
causes a pure rotation of p about an arbitrary axis by an arbitrary angle 
(3 degrees of freedom), a "poldistortion", which exchanges an arbi- 
trary amount of power between Stokes I and and an arbitrarily oriented 
component of p (3 degrees of freedom), and an overall system gain 
(one degree of freedom). A similar decomposition was made indepen- 
dently by Britton 1 2000), who notes that Jones matrix transformations 
are isomorphic with the Lorentz group, wi th the type of transforma- 
tion called a "poldistortion" bv iHamaked fcOOd) corresponding to a 
Lorentz boost. Apart from providing a convenient geometric interpre- 
tation of the kinds of transf ormations effected by linear system com- 
ponents, as Hamaker (2000) notes this decomposition suggests a cali- 
bration strategy that begins with the observation of a presumed unpo- 
larized source, immediately giving the gain and poldistortion. The pol- 
rotation can then be obtained by observations of a further two sources 
about which some assumptions can safely be made (for example, that 
they are linearly polarized and/or have a given position angle). 
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Unfortunately efforts to calibrate the WSRT in tied array mode 
are severely hampered by the inclusion of non-linear system compo- 
nents. Specifically, the voltage signals from each telescope are sam- 
pled using thresholds that are dynamically determined using a control 
system that attempts to maintain the variances of the sampled signals 
at a constant level, based on averaging with some time constant of 
the order 1 ms. These are known as the automatic gain controllers 
(AGCs). Although these components might at first appear disastrous 
for the detection of non-stationary signals, the fact that they normalise 
the individual telescope powers while the source signal adds coher- 
ently means that the distortion is reduced for N » l 3 . However, the 
fact remains that the AGCs introduce an unknown, time-varying gain 
to each of the real and imaginary parts of each polarization channel 
of each telescope, making precise calibration impossible for tied array 
observations, and we strongly suggest that future telescope arrays pro- 
vide tied-array systems free of AGCs. Since the system is not precisely 
calibratable, we limit the scope of our calibration to the grossest of in- 
strumental effects, differential gain and phase between the summed X 
and Y polarization channels. Remaining effects are expected in any 
case to be small, d ue to the high accuracy of dipole setting at WSRT 
IWeiler & Wilsonlll977l) . and the tendency for random error terms to 
cancel in the sum of Eq. IA.4l By comparison of our observations with 
published polarimetric profiles available via the EPN database 4 , we 
estimate that our results are accurate to within a few percent, certainly 
sufficient for examination of the basic polarization properties of the 
source. 

The procedure we use for determining the differential phase and 
gain is based on the differential Faraday rotation across the observ- 
ing band. This provides a variety of input polarization states incident 
on the telescope, with a known relationship between them, allowing 
one to fit simultaneously for the source polarization and the instru- 
mental response. In this regard the technique is similar to previous 
strategies employing the parallactic angle rotation during long obser- 
vations of a source IStinebring et alll984al) . which cannot be used at 
the WSRT because its equatorial mounts cause the dipoles to track 
the parallactic rotation (by the same token, this eliminates a potential 
source of time- variability in incompletely calibrated measurements). 
Just as the parallactic technique is limited to the assumption that nei- 
ther the source polarization nor the telescope response changes with 
time, our method assumes frequency independence except for an over- 
all factor incorporating the gain and intrinsic intensity spectrum of the 
source. The technique is also limited by the commutativity of certain 
transformations with the known transformation effected by Faraday or 
parallactic rotation, about Stokes V. That is, for any solution of the 
system Mue ller matrix (whi ch is constructed from the Jones matrix, 
see e.g. Hamaker et al 119961) satisfying 



Si = MR v (29i)S, 



(A.5) 



where S is the intrinsic Stokes 4-vector, S,-' is a measured (uncali- 
brated) Stokes 4-vector observed under a known Faraday or parallactic 
rotation angle and R v (20j) is a Mueller matrix effecting a rotation 
of angle 26>, about Stokes V, there exists an infinite number of other 
solutions 



= MER v (28j)E~ [ S, 



(A.6) 



where E is any matrix that commutes with Ry{26j). Without addi- 
tional constraints, this approach is therefore only able to determine 
the system response to within an arbitrary rotation ("polrotation") and 



3 See the PuMa manual, 

http : / /www . astron . nl/wsrt/WSRTobs/PuManual . pdf 

4 http : //www . mpifr-bonn . mpg . de/di v/pulsar/data/ 
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Fig. A.l. Observed values of Stokes / (open circles), Q (trian- 
gles), U (squares) and V (filled circles) in the peak longitude 
bin of PSR B0329+54, as a function of frequency. Predictions 
of the best-fitting model (incorporating all longitude bins) are 
shown with lines. 



boost ("poldistortion") in Stokes V. However, as stated above we have 
chosen to use a reduced description of the system response, parame- 
terised only by frequency-dependent gain, and X-Y differential gain 
and phase. For assumed perfect linear dipoles, differential phase and 
gain errors c orrespond to rotations and boosts about Stokes Q (e.g. 
Britton 2000), which do not commute with R y, obviating the need for 
additional constraints. Our model is thus: 



S'(v) = g(v)R Q (2RM(c/v) 2 + e)B Q (J3)S + noise, 



(A.7) 



where g(v) models the frequency dependent gain and source spec- 
trum, RM is the rotation measure of the propagation path between 
the source and Earth, 8 is the phase offset between X and Y channels, 
and p = Ingx/gy where g x and g y are the (voltage) gains of the X and 
Y channels. 

We define the best-fit system response as the one that minimises 
the global^ 2 : 



z 



(£1 A| 

2 2 1 2 

I,i Q,i U,l V,l 



(A.8) 



where i indexes all measured Stokes vectors (i.e. all frequency chan- 
nel/pulse longitude pairs), cr^jQ UVj are the variances of each of the 
Stokes parameters in off-pulse regions (potentially also dependent on 
z), and 



Di=S'i- MiRyWdS, 



(A.9) 



where in our case, Mi = g(Vi)R e (2RM(c/v,) 2 + 6)B Q (P). In the gen- 
eral case where a full system calibration is achievable M should be 
parameterised as a function of frequency by the seven real and imag- 
inary components of the elements of the Jones matrix (the imaginary 
part of one of the diagonals can be set to zero since absolute phase 
is irrelevant). Johnston (2002) describes a similar scheme to that de- 
scribed above, employing a Jones matrix formalism and parallactic 
angle rotation, from which we differ in several respects. Firstly rather 
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than minimising x 2 , which corresponds to the maximum likelihood 
estimator, Johnston (2002) minimises the Manhattan distance, which 
is likely to have poorer esti mation accuracy than maximum likelihood. 
Secondly, Johnston (2002) discusses the difficulties of multi-step cal- 
ibration schemes but does not propose an alternative. Including all 
measurements of all sources in the global^ 2 , including any injected 
noise sources and continuum calibrators, and minimising it in one 
step is in our view the best way to avoid these difficulties. Finally, 
we caution against neglecting second-order terms in the expansion of 
the Mueller matrix as it is unnecessary and introduces the possibility 
for error if and when the technique is applied to telescopes with non- 
negligible cross terms (e.g. LOFAR/SKA). We note that the issues de- 
sc ribed above have re ceived independent treatment in the recent work 
of Ivan Stratenl l2004) . who also provides a method for minimization 
of Eq. lA,8l that avoids explicit expansion of the Mueller matrix and its 
partial derivatives. 

We performed the procedure described above using a frequency- 
resolved polarimetric profile from the first 12500 pulses of the ob- 
servation of PSR B0329+54 described in the text. That 6 =s ±90° is 
clearly visible in the frequency dependence of the Stokes parameters 
(Fig. IA.lt . where the Faraday modulation appears mainly in Q and V, 
instead of Q and U as expected if 6 = 0. In fact it was this feature, seen 
in this and other WSRT observations in the 328 MHz band that led us 
to the calibration scheme described here. The results of the fit support 
the assertion that 9 t 0, yielding 6 = -75.24° and fi = -0.0278, im- 
plying gx/gy = 0.973 (the model predictions using these values and the 
other parameters are plotted in Fig. lA.U . These values were used to 
correct the recorded Stokes vectors before further analysis described 
in the body of this paper. 



